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ABSTRACT 

We compnte a series of three-dimensional general relativistic magnetohydro¬ 
dynamic simnlations of accretion flows in the Kerr metric to investigate the prop¬ 
erties of the nnbonnd ontflows that resnlt. The overall strength of these ontflows 
increases sharply with increasing black hole rotation rate, bnt a nnmber of generic 
featnres are fonnd in all cases. The mass in the ontflow is concentrated in a hollow 
cone whose opening angle is largely determined by the effective potential for mat¬ 
ter orbiting with angnlar momentnm comparable to that of the innermost stable 
circnlar orbit. The dominant force accelerating the matter ontward comes from 
the pressnre of the accretion disk’s corona. The principal element that shapes 
the ontflow is therefore the centrifngal barrier preventing accreting matter from 
coming close to the rotation axis. Inside the centrifngal barrier, the cone con¬ 
tains very little matter and is dominated by electromagnetic helds. The magnetic 
heldlines inside the cone rotate at a rate tied closely to the rotation of the black 
hole, even when the black hole spins in a sense opposite to the rotation of the ac¬ 
cretion flow. These helds carry an ontward-going Poynting hnx whose immediate 
energy sonrce is the rotating spacetime of the Kerr black hole. When the spin 
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parameter a/M of the black hole exceeds ~ 0.9, the energy carried to inhnity by 
these outflows can be comparable to the nominal radiative efficiency predicted 
in the Novikov-Thorne model. Similarly, the expelled angular momentum can be 
comparable to that accreted by the black hole. Both the inner electromagnetic 
part and the outer matter part can contribute in signihcant fashion to the energy 
and angular momentum of the outflow. 

Subject headings: Black holes - magnetohydrodynamics - jets - stars:accretion 


1. Introduction 

Whether from supermassive objects in galactic nuclei or stellar mass objects in Galactic 
binaries, outflows—often relativistic—are frequently seen from black holes. However, there 
is little that is known with conhdence about their nature and dynamics. The essential 
ingredients are, of course, well-known: magnetic field, accretion, and rotation. But basic 
questions such as the forces that drive the jets, the mechanisms that regulate their content, 
and the constraints that collimate them remain open. Magneto-centrifugal effects have 
received the greatest attention (Blandford & Znajek 1977; Blandford & Payne 1981; Shibata 
& Uchida 1985; Punsly & Coroniti 1990), but there are many possible variations on this idea. 
In most such studies, magnetic held boundary conditions are guessed and their consequences 
derived without much consideration of how the posited held structures would arise in the 
context of accretion hows. 

We approach these questions from the opposite point of view: we simulate numerically 
how accretion dynamics self-consistently create magnetic held and explore the outhows that 
result. Employing the general relativistic three dimensional MHD code described in De 
Villiers & Hawley (2003), we have previously investigated the general character of accretion 
onto black holes in a series of simulations for which an overview is presented in De Villiers, 
Hawley & Krolik (2003; hereafter Paper I). The magnetic held structures that develop in 
these accretion hows were described in Hirose et ah (2004; hereafter Paper H). We discussed 
the resulting outhows in De Villiers et al. (2005; hereafter Paper HI), and the dynamics 
of the inner disk in Krolik, Hawley, & Hirose (2005; hereafter Paper IV). In this paper we 
return to a more in-depth consideration of the type of unbound outhows that were discussed 
in Paper HI. We will present a more quantitative and dynamically-oriented analysis of the 
outhows, and draw on several new simulations designed to elucidate their properties. We 
will focus our attention on three issues in particular: the underlying principles governing 
the material component in the outhow; the dynamics and spin-dependence of the outhow’s 
electromagnetic segment; and the relationship between the electromagnetic contribution and 



the class of ideas placed under the heading of the Blandford-Znajek mechanism. 


2. Overview of Simulations 
2.1. Qualitative character 

In the present work we solve the equations of ideal MHD in the metric of a rotating 
black hole, while assuming that the gas is adiabatic except in regions of strong compression, 
where an artihcial bulk viscosity is employed to capture the entropy produced in shocks. 
The specihc form of the equations we solve, and the numerical algorithms incorporated into 
the GRMHD code, are described in detail in De Villiers & Hawley (2003). 

In this paper, we present several new simulations, which, like the simulations reported 
in Papers I-IV, follow the evolution of an isolated magnetized gas torus orbiting near a 
black hole. The new simulations differ from the previous ones in part because of several 
small technical changes in the code that permitted generation of new, improved versions of 
previously-published simulations (e.g., raising the cap on the Lorentz factor, saving addi¬ 
tional data, extending the radial grid closer to the horizon for a/M = 0.9). In addition, we 
have studied several new high-spin cases: a/M = 0.93, 0.95, and 0.99, and a case with a 
counter-rotating black hole, a/M = —0.9. 

The late-time structure of all our simulations contains the same basic elements pointed 
out in Paper I, although their relative importance and location can change considerably as 
a function of black hole spin. These structures are illustrated in Figure 1, a schematic of the 
accretion system as outlined in Papers I through IV. Along the equator there is a somewhat 
thickened, wedge-shaped Keplerian accretion disk. A net accretion flow is produced within 
this disk by MHD turbulence whose ultimate origin is the magneto-rotational instability 
(Balbus & Hawley 1998). In the disk, the magnetic held is subthermal; the ratio of the gas 
to magnetic pressures, i.e., Pgas/Pmag = P, is on average on order 10 to 100. A short distance 
outside the marginally stable orbit around the black hole, the equatorial pressure and density 
reach a local maximum in a region we refer to as the inner torus. The location and vertical 
thickness of this inner torus can and does vary as a function of time throughout a simulation. 
Inside this local pressure maximum, the density and pressure drop as the how accelerates 
inward toward the black hole in the plunging region. The plunging region typically begins 
near the innermost marginally stable circular orbit of the black hole, The high latitude 
region above and below the disk is the corona, a region of hot, magnetized plasma with a 
magnetic held whose strength is, on average, near equipartition (/5 ~ 1). Along the spin axis 
of the black hole there is a centrifugal funnel that is largely empty of matter, but hlled with 



4 


magnetic field. This magnetic field transports energy outward in the form of a Poynting 
flux. There is a region of unbound mass flux at the boundary between the evacuated funnel 
and the corona which we refer to as the funnel wall jet. This funnel and its boundary are, 
of course, the center of attention in this paper. 


2.2. Numerical technique 

The GRMHD code we employ for this study solves the equations of ideal MHD in 
the metric of a rotating Kerr black hole of mass M in three spatial dimensions expressed 
in Boyer-Lindquist coordinates, (t,r,9,(l)). We use the metric signature (—,+,+,+). The 
determinant of the 4-metric is g, = a ^ 7 , where the lapse function a = and 

7 is the determinant of the spatial 3-metric. We follow the usual convention of using Greek 
characters to denote full space-time indices and Roman characters for purely spatial indices. 
We use geometrodynamic units where G = c = 1; time and distance are in units of the black 
hole mass, M. The spin of the black hole is specihed by a, with a/M < 1. 

The primitive variables that describe a relativistic test fluid are its density p, specihc 
internal energy e, four-velocity U^, isotropic pressure P, and magnetic held four-vector b^. 
The relativistic enthalpy is h = l + e + P/p. The pressure is related to p and e through 
the equation of state of an ideal gas, P = pe (T — 1), where T is the adiabatic exponent. 
For these simulations we take T = 5/3. The code evolves auxiliary variables constructed 
from the relativistic boost factor W = namely D = pW, E = pehF, a four-momentum 
= {ph + ll&ll^) hFand a transport velocity W = W/U*. The magnetic held of the 
huid is evolved using the constrained transport magnetic held, Fjk = [ijk] B\ from which 
the magnetic held four-vector is constructed. The ideal MHD condition 

requires = 0. 

The algorithms employed in this GRMHD code resemble those of the well-known ZEUS 
code (Stone & Norman 1994a,b) used for nonrelativistic astrophysical MHD simulations. 
The evolution is divided into a transport step, where quantities are advected through the 
grid by the velocities U*, and a source step where forces are applied. We use the Gonstrained 
Transport approach for the induction equation (Evans & Hawley 1988) so that the magnetic 
held remains divergence-free. The evolution timestep is set by the smallest light-crossing 
time through a grid zone. 

In contrast to some alternative numerical schemes that have been developed for general 
relativistic MHD (e.g., Gammie, McKinney & Toth 2003; Komissarov 2004) our scheme is 
not fully conservative. Rather than evolve the stress-energy tensor directly, we work with 
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Fig. 1.— Schematic illustration of the main dynamical features seen in the accretion disk 
simulations. The jet is dehned as those regions of unbound positive radial momentum. 
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the auxiliary variables from which the primitive variables are easily recovered. We evolve 
only the three spatial components of the four-momentum and use velocity normalization to 
obtain the St component. The boost factor W is then computed using but because the 
dehnition of the inertia used in the code includes W, it is found from the solution to a quartic 
equation derived from four-normalization (De Villiers, Staff, & Ouyed 2005). The internal 
energy is evolved independently. An artificial viscosity is included to thermalize some of the 
kinetic energy lost in shocks, but there is no mechanism included to capture magnetic energy 
lost through magnetic reconnection. 

It is worthwhile to review briefly some of the known strengths and weaknesses of the 
GRMHD code. In tests and applications carried out to date the code has proven to be 
fast and robust, capable of carrying out long-term simulations in three spatial dimensions. 
Indeed, as of the time of the writing of this paper, this is the only code that has been able 
to perform long-term (several thousand M in time) three-dimensional GRMHD simulations. 
The code handles both the weak held (i.e., high f3 = Pgas/T’mag) limit within the accretion 
disk and the strong helds {j3 1) found in the evacuated funnel. Gomparisons between 

simulations using nonrotating black holes and nonrelativistic simulations using the ZEUS 
code show good qualitative agreement within the region where gravity is well approximated 
by a pseudo-Newtonian potential. Our simulations also show qualitative agreement with 
the axisymmetric simulations of black hole accretion carried out by McKinney & Gammie 
(2004). 

One of the limitations of the code is that there must be a minimum hoor on the huid 
density, set in the simulations of this paper at 10“^° of the initial density maximum. This 
density hoor value is chosen simply to provide a wide dynamic range between the density 
within the matter-dominated regions and the region that should be regarded as the numerical 
“vacuum.” To prevent possible numerical instability in near-vacuum zones we place a cap 
on the value of the boost factor; here Wmax = 10 or 6. Because we evolve an internal energy 
equation, the thermodynamics are only approximate. Some amount of energy that would 
otherwise be thermalized (e.g., in magnetic reconnection) is lost. As we have no explicit 
cooling, the disks grow hot as they accrete, but not as hot as they would be with complete 
energy conservation. Gonversely, of course, they are not as cool as they might be if radiation 
losses were taken into account. 

Because we use Boyer-Lindquist coordinates, the horizon represents a coordinate singu¬ 
larity, and our inner grid boundary, must he outside this point. The metric terms increase 
rapidly near the inner boundary, making accurate numerical calculations there more difficult. 
In practice we select an inner radius such that the Kerr metric term A = — 2Mr + {a/M)‘^ 

is approximately ~ 0.1-0.02. As a result, the spacetime rotation frequency at rin is less than 
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that of the horizon. The inner boundary uses zero gradient conditions on the fundamental 
variables, and the radial velocity and momentum at are not allowed to become positive. 

Even apart from the coordinate singularities at the horizon in Boyer-Lindquist coordi¬ 
nates, it is well known that spherical coordinates have a singularity along the axis. As a 
practical matter, the closer one is to the 0 = 0 or tt axis, the more narrow and wedge-like the 
(j) zones become, leading to a substantial reduction in the Courant-limited timestep among 
other problems. To avoid this, we set the 6 boundary not at the axis but at a small angle 
off the axis and employ reflecting boundary conditions. This means that with our grid we 
cannot follow a flow passing over the north or south pole of the black hole. 

These conditions imposed on the computational grid will, of course, affect the results. 
The “torque” applied by the rotating spacetime will be limited to that provided exterior 
to the inner boundary. Tin, and there will be no jet power along the cutout on the axis. 
The numbers reported for the models should not be regarded as quantitatively predictive of 
the jet power a black hole with a specihc a/M value would have in Nature. Rather we are 
interested in the relative values of quantities across the full complement of the simulations 
and what that can mean for the role of a rotating hole in the combined accretion plus jet 
production process. 

It is important to keep these code limitations in perspective. As no numerical scheme is 
perfect, all simulators must make choices. There are often alternatives that can address one 
specihc limitation, but introduce new limitations of their own. For example, the problems 
associated with a spherical mesh can be eliminated through the use of a Cartesian grid, but 
then angular momentum is no longer a directly conserved fundamental variable and main¬ 
taining a simple orbit becomes a matter of carefully balancing forces acting on distinct linear 
momentum components. Diffusion error associated with motion across the grid can then be 
a signihcant source of nonphysical viscosity in an accretion disk. When evolving a total 
energy equation, the internal energy is dehned as what remains after all other independently 
evolved energies (kinetic, magnetic) have been subtracted off. If the internal energy is a small 
fraction of the total energy (as it is whenever h/r 1), its evolution can be dominated by 
truncation errors in the kinetic and magnetic components. In relativity, the derivation of 
primitive variables from conserved quantities is difficult and time-consuming. 

At this stage the only practical approach is to simulate well-chosen problems with a 
variety of numerical codes and to evaluate the outcomes both in comparison to other results 
and with self-consistency checks, keeping a specihc code’s limitations hrmly in mind. 



2.3. Problem specifications 


For this paper we have computed a new set of models in a spherical quadrant compu¬ 
tational domain. The inner radial boundary, denoted rjn, is placed just outside the event 
horizon (the exact amount varying from case to case), and the outer radial boundary is at 
^out = 120M. The azimuthal angle ranges over 0 < 0 < 7r/2 and the polar angle over 
O.OdSTT < 6 < O.QSStt. The radial boundary conditions are pure outflow, the azimuthal 
boundary conditions are periodic, and the boundary conditions near the polar axis are re¬ 
flecting. The grid consists of 192 radial cells spaced using a hyperbolic cosine function in 
order to maximize the resolution near the inner boundary, 192 polar angle cells, with an 
exponential grid spacing function that concentrates zones near the equator, and 64 evenly- 
spaced cells in azimuthal angle. This resolution is close to the maximum feasible at this time 
for computing an ensemble of three dimensional simulations. 

These models are of the same form as the KD models described in Paper I. Since those 
simnlations were run, we have added new run-time diagnostics and have stndied a greater 
range of black hole spin parameters: a/M = —0.9, 0., 0.5, 0.9, 0.93, 0,95, and 0.99. Note 
that we have included one case with a black hole that is connter-rotating with respect to 
the disk. As in Paper I, the initial conditions consist of an orbiting isolated gas torus with 
a pressure maximum located at r ^ 25M, and a slightly snb-Keplerian initial distribntion 
of angnlar momentum throughout. We construct an initial torus so that its inner edge is 
located at r = 15M regardless of the black hole spin rate. The initial magnetic held consists 
of loops of weak poloidal held lying along isodensity snrfaces within the torns. The strength 
of the magnetic held is chosen so that the ratio of the total magnetic pressure to the initial 
total gas pressnre in the torus is 0.01 {(3 = 100). Table 1 lists the models along with some 
reference valnes such as the radius of the marginally stable orbit, r^s, and the frame-dragging 
rate at the inner grid bonndary, a;(rin). 


3. Global Measures of the Outflow 

The strength of the outflow may be characterized by how much mass, linear momentum, 
angular momentum, and energy it carries to large distance. We compnte these measnres 
throngh shell-integrated flnxes that are calculated using 

{J^){r,t) = J J Q{r)^/^ded(j), (1) 

for qnantity Q at radius r, where the bounds of integration range over the full 6 and 0 com¬ 
putational domains. These fluxes are computed and stored every M in time, and integrating 
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the flux over time yields the total outflow or accretion rate. Whenever we speak of mean 
properties of the outflow, we require a criterion to distinguish those cells participating in 
the outflow from those in the accretion flow. Our definition of “outflow” is any cell with 
— hUt > 1 (i.e., that is unbound) and with f/” > 0 (i.e., moving outward). This criterion 
works well for these simulatons, as only the outflow near the axis (the “jet” outflow) is 
unbound; the backflow from the disk itself remains bound. In practice we divide the shell 
integrals into two parts, one for bound and one for unbound flow. We will indicate the shell 
integrals specific to either the outflow or the accretion flow by multiplying the integrand by 
either W or 1 — W, where W is a function whose value is unity in outflow cells and zero 
otherwise. 

In the present work we are concerned with the strength of the jet compared to what is 
accreted into the hole. We would like to characterize the outflows in terms of the accretion 
rate, but this can be ambiguous because of the scale-free nature of these simulations (code 
values of density are arbitrary), and the fact that the total mass of the initial torus is a 
function of black hole spin. We proceed by computing the total mass, angular momentum 
and energy (in code units) that leave the grid going into the black hole; these are the “hole” 
values. Second we compute the total unbound mass, angular momentum and energy that pass 
through r = lOOM; these are the “jet” values. To focus on the quasi-steady state portion of 
the simulation, we compute these totals from t = 3000 to 10‘^M. These integrated quantities 
are listed in Tables 2, 3 and 4, along with measures of the outflow efficiency. For the case 
of the angular momentum and energy fluxes, we further distinguish between that carried by 
the matter (superscript m) and that carried by the electromagnetic fields (superscript em). 
The efficiencies are defined as 


= (E™ - Mjet) /Mhole (2) 

^ F^fe?/Mhole, (3) 

where Mhoie is the total rest-mass accreted into the black hole, while Mjet and Ej^t are 
the rest-mass and energy carried off through the r = lOOM surface. Rest-mass energy is 
subtracted from the matter portion so that the efficiency refers only to “usable” energy. In 
addition to relative efficiency, it is also useful to compared to the standard Novikov-Thorne 
(Novikov & Thorne 1973) thin disk efficiency, which is the binding energy of the innermost 
stable circular orbit. These efficiencies are also listed in Table 4. 

Before discussing these numbers, we must first express a caution: the fluxes are, in 
general, not constant in time or radius. Total (i.e., electromagnetic plus matter) quantities 
can depend on radius, even in the long-term average. Even where the total is constant, 
electromagnetic fields and matter can, and do, exchange energy and momentum. In most 
cases, for example, the EM outflow rates reach a peak at several tens of gravitational radii 
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and decline slowly from there to the outer boundary. In contrast, the magnitude of the time- 
averaged rest-mass expulsion rate increases monotonically toward the outside, although its 
rate of increase decreases with growing distance from the center. In a strict sense, the outflow 
numbers apply only at our hducial outflow radius r = lOOM. Consequently, the numbers 
we give in the tables should be regarded as useful for comparing relative changes between 
simulations, and to study how outflow properties depend on black hole spin, not as providing 
some sort of absolute measure of jet power. 

Examining the table of mass and energy flows, the most obvious conclusion is that for 
rapidly-rotating black holes the jet efflux is not insignificant compared to what is accreted. 
In other words, the jet plays a noticeable role in the global energy budget of accreting black 
holes. The mass outflow is all within the funnel wall jet, whereas the Poynting flux outflow is 
primarily within the funnel itself. The ratio of jet mass to accretion mass is generally a few 
percent, although it rises steeply from a/M = 0.95 to a/M = 0.99. Moreover, although the 
trend in these simulation data is not monotonic, there is a general tendency for the expelled 
energy to increase sharply with increasing black hole spin. The outflows are almost nil in 
the Schwarzschild case, weak when a/M = 0.5, and become substantial when \a/M\ > 0.9. 

When \a/M\ > 0.9, the total jet efficiency is generally several tens of percent, with the 
matter portion somewhat larger than the electromagnetic part. The retrograde case is less 
efficient in driving a jet than the prograde case with the same magnitude of spin, but only by 
a factor of two. Interestingly, in otherwise similar simulations most of whose results closely 
resemble ours, McKinney (2005) remarks that he finds the retrograde jet efficiency to be 
smaller by a factor of ten. 

The ratio of angular momentum in the unbound outflow to that deposited in the black 
hole is also a very strong function of black hole spin. Only a few tenths of a percent in the 
Schwarzschild case, it rises to ~ 25% when a/M = 0.9. In the highest spin case reported here 
{a/M = 0.99), the net angular momentum accreted is extremely small because the rate at 
which the black hole spins down via electromagnetic torques nearly matches the rate at which 
it acquires angular momentum by accreting matter. Consequently, the angular momentum 
expelled in the outflow is an order of magnitude greater than the net amount captured. In all 
cases, the electromagnetic portion of the angular momentum carried outward is comparable 
to the matter portion carried in the funnel wall jet. 

Note in particular that the electromagnetic angular momentum carried outward in the 
retrograde case has the opposite sign relative to the accreting matter’s orbital angular mo¬ 
mentum. One consequence of this is that the net angular momentum outflow in this case is 
very small because the electromagnetic part very nearly cancels the matter part. This fact is 
also an indication (that will be confirmed in the following section) that the electromagnetic 



outflow is directly tied to the black hole’s rotation. 


4. Electromagnetic Outflow: the Funnel Interior 
4.1. Character of magnetic held in the funnel 

In the initial state there is no magnetic held at all within what eventually becomes the 
outflow region. All of the simulations presented here rapidly develop a large scale held within 
the funnel, which in turn leads to the formation of the Poynting flux jet. This is one of the 
signiflcant results of simulations of this type (see also Paper I and II, as well as McKinney 
& Gammie 2004), namely that a large scale, jet-driving magnetic held can develop naturally 
as the result of magnetically controlled accretion. 

This funnel held formation process may conveniently be visualized in terms of three 
steps. These are illustrated in the panels of Fig. 2, which show the evolution of the (3 
parameter in the inner region of the inflow in run KDPg. First, as the magnetorotational 
instability (MRI) develops in the disk, matter begins to move inward. The magnetic held is 
stretched radially, and this radial field is sheared into an increasingly strong toroidal field. An 
extremely thin precursor of matter moves quickly inward in the equatorial plane, surrounded 
by a thicker envelope of magnetically-dominated {j3 ~ 0.1) matter. The matter density above 
and below the equatorial accretion flow is very low, in part because of the initial conditions. 
However, the expectation of low density above and below the plunging inflow is likely to be 
general because it is the centrifugally-excluded zone. When the accretion flow reaches the 
black hole, matter can drain off the field lines. The vertical component of is then great 
enough to drive a vertical expansion of the held (center panel in Fig. 2). A “magnetic tower” 
forms, as has previously been predicted on analytic grounds (Lynden-Bell 2003) and seen in 
other simulations (Kato, Mineshige, & Shibata 2004). It lasts for a time ~ 100-200M. 

The steady-state structure is created in the third step. Because the top of the “magnetic 
tower” expands more rapidly than its base, the toroidal wrapping of held lines stretches into 
a less tightly-wound helix. To stretch the held from the “tower” region ~ 20M from the 
plane all the way to the outer boundary requires only ~ lOOM in time because the motion is 
relativistic. When this happens in a Schwarzschild metric, the held lines end up very nearly 
radial; with increasing black hole spin, the field lines wind more tightly (Paper II). That is, 
it is only the intrinsically general relativistic effect of frame-dragging that in the steady state 
imparts a net helical twist to the field far from the black hole. Any remnant of the “tower” 
that persists in the long-run is not due to MHD accretion dynamics; it is, rather, an imprint 
of the rotating metric. 
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Fig. 2.— Ratio of the azimuthally-averaged gas to magnetic pressure (3 aX t = 560M, 640M, 
and 720M in KDPg. The color contours are in a logarithmic scale: dark red is gas pressure 
dominated (/9 = 10^) and dark blue is magnetic held dominated (/9 = 0.001). 
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As discussed in Paper II, the interior of the outflow funnel is magnetically-dominated 
in the sense that \\h\\^/{ph) 3> 1. Values of this ratio in the range ~ 10-10^ are common. 
However, we stress that the particular value of ph found in the funnel is not well-determined 
in these simulations. The density reaches the floor value at late time in much of the funnel. 
This occurs for all values of the floor that have been tested; clearly the funnel is being 
evacuated. Except for some low density initial background material, the matter on the grid 
all has angular momentum, and this prevents it from entering the axial funnel. Thus the 
Ending that the funnel should be evacuated is expected, and our results should simply be 
taken as saying that the magnetic energy density is likely to be at least several orders of 
magnitude greater than ph. 

The steady-state funnel held is illustrated in Figure 3, which shows the ratio of radial 
to azimuthal components of the magnetic held in the outhow at a late time in three of the 
simulations. Because the gttp term in the metric becomes small when r 3> M, and the focus 
here is on behavior at large distances, we have approximated this ratio as if the metric were 
purely diagonal, i.e., we have plotted the ratio /B^) {grr!■ As the hgure shows, the 
held is very purely radial when the black hole does not spin, but grows progressively more 
toroidal with increasing black hole rotation rate. For intermediate values of spin, the held is 
more nearly radial near the polar axis, more toroidal at larger angles. 
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Fig. 3.— Azimuthal mean of the magnetic radial to azimuthal component ratio at a late 
time in the a/M = 0, a/M = 0.5, and a/M = 0.9 simulations, from left to right. The 
contours are logarithmic. 
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If this flow were time-steady and axi-symmetric, it would make sense to discuss the rate 
at which the poloidal portion of the held rotates around the axis. Under those conditions the 
heldline rotation rate is dehned as cn = or equivalently as cn = — B'^V^/B^. 

Here, because the situation is neither axi-symmetric nor stationary, the concept of heldline 
rotation rate is not so clear. Nonetheless, in order to provide some intuition by comparison 
to a more familiar context, we construct a local and instantaneous dehnition of heldline 
rotation rate: 


UJ = 


_ , V^B^9rr + V^B^gee 
{B^Ygrr + m^gee' 


( 4 ) 


The second term describes how fast matter slips along the toroidal held relative to the 
rotation of the poloidal part. In the time-steady axi-symmetric limit, equation 4 reduces to 
the standard dehnition of heldline rotation (as found, e.g., in Phinney 1983). 


In Figure 4 we show how this dehnition of the rotation rate compares to specihc predic¬ 
tions made in the context of the Blandford-Znajek model. When, for example, the heldlines 
are radial, it is expected that the rotation rate close to the polar axis should be half the black 
hole rotation rate (Blandford & Znajek 1977; Phinney 1983). To make this comparison, we 
must average over both position along individual heldlines and over time. For the former 
average, we make use of the fact that the how is nearly conical at large radius, and dehne 


{uj){t,e) 



dru(t, r, 9, 0), 


( 5 ) 


where we choose rjn = rout/2 in order to restrict the average to the nearly conical region. In 
practise, although we see substantial huctuations in cn as a function of radius at hxed 9, we 
see little in the way of systematic trends with radius in this outer region, so the exact range 
of averaging should not make much diherence. 

The details diher slightly as a function of time within any one simulation and also from 
simulation to simulation. For example, the behavior of the held line rotation at the axis 
boundary varies, sometimes dropping oh as in Figure 4 and sometimes increasing. (Recall 
that the axis has a rehecting boundary condition applied a.t 9 = 0.0457r.) Within the main 
part of the funnel, however, things are more consistent from model to model. The heldline 
rotation rate is close to, but generally ~ 20-25% smaller than, half the spacetime rotation 
rate of the inner boundary (see Table 1). We use this criterion rather than half the rotation 
rate of the black hole itself because the latter quantity plays no direct role in our simulation. 
The diherence between the classical prediction and the observed rotation rate in our results 
is larger than that observed in the two-dimensional simulations of McKinney & Gammie 
(2004). They found that, within the force-free part of the outhow, the rotation rate was 
about 10% less than half the black hole rotation rate. 




e 

Fig. 4.— Fieldline rotation (as defined in text) in the outfiow for the simulation with 
a/M = 0.95. The horizontal lines show half the rotation rate of the black hole (dotted) and 
half the rotation rate of the inner boundary (dash-dot). 
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We further find that the rotation rate of the funnel heldlines is always in the black hole’s 
sense; that is, in the retrograde case, the heldlines rotate in the same direction as the black 
hole, but in the opposite direction of the accreting matter. Specihcally, dehning the black 
hole spin as in the positive direction, when the disk and the black hole rotate in opposite 
senses, the rotation of the funnel-wall matter is backward, and, at large radius, of smaller 
magnitude than the heldline rotation. This contrast nnderlines the fact that the heldlines’ 
rotation is fundamentally a result of frame-dragging in the ergosphere. Interestingly, the 
mean rotation rate in the retrograde case is actually considerably closer to half the rotation 
at the inner boundary than it is in the prograde cases. 

Finally, althongh we have hopes of hnding large-scale held-generation processes that are 
completely general, we note that the fnnnel held formation seen here does depend in certain 
ways upon the problem conhguration. First, the crncial ingredient is the presence of radial 
held extending throngh the plnnging region as material hrst accretes into the hole. The 
polarity of this held, and of the snbseqnent fnnnel held, is a conseqnence of the initial held 
conhgnration, namely poloidal loops with one sign of the radial held above the eqnator and 
the opposite sign below. A simnlation with a purely toroidal initial held (discussed in Paper 
III) is a connter-example that does not form an organized funnel held. Second, the initial 
plnnging inhow is snrrounded by a near-vacunm which ohers no opposition to the formation 
of the magnetic tower. Third, the accretion how is relatively hot, with H/R ~ 0.15. It is 
possible that thinner, colder disks would not generate sufficiently strong helds to create the 
magnetic tower; additional simulations will be necessary in order to explore the degree of 
generality of this large-scale held-generation mechanism. 


4.2. Angular distribution of Poynting flux 

In the time-steady, axi-symmetric limit in which it is possible to dehne clearly a heldline 
rotation rate, the Poynting hux on the inner bonndary can be described in terms of the 
magnetic held strength and the heldline rotation rate: 


R = sin^ 9, (6) 

where H = R + (a/M)^ cos^0 and A = R — 2Mr + {a/MY are the standard fnnctions 
appearing in the Boyer-Lindqnist metric (Blandford & Znajek 1977). Although A —0 
as the event horizon is approached, the magnetic held strength diverges in an appropriate 
way to make T continuous. In the specihc circumstances envisaged by Blandford & Znajek 
(1977), the polar angle dependence of the magnetic held intensity follows from the condition 
of satisfying a force-free eqnilibrinm for a held of specihed multipolar character. They 
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predict that for a split monopolar field it should diminish by about 20% from the axis to 
the midplane. In addition, they predict that the heldline rotation rate should be close to 
1/2 the black hole rotation rate almost independent of polar angle, so the overall polar angle 
dependence of the Poynting flux should be almost oc sin^0. 

In our simulation, we hnd that the magnetic held sufficiently dominates the matter 
inside the funnel that a force-free conhguration might be expected on average (although 
certainly not obtained instantaneously, as the huctuations are large), while the detailed 
dynamics of the coupled accretion and jet hows determine the geometry of the held. In the 
simulations with \a/M\ < 0.9, the time-averaged distribution of T (calculated directly from 
the electromagnetic contributions to T/, not from equation 6) does indeed vary with polar 
angle very nearly as sin^ 0 (Fig. 5). However, this equilibrium structure is not necessarily 
realized at any given time, for there are tens of percent huctuations around the mean value. 
In the larger spin cases, these huctuations are larger in magnitude. Even after averaging over 
some 27 snapshots, the polar angle dependence of the Poynting hux in the a/M = 0.95 case 
is still quite rough. It is also noteworthy that, although the Poynting hux goes inward in the 
midplane accretion how at lower spins, for a/M > 0.95, it is outward almost everywhere, 
even within the equatorial accretion how. The overall pattern seen in the a/M = 0.95 case— 
a rise from the polar axis to mid-latitudes, a sharp drop at an angle ~ 7r/4 from the pole, 
and irregularly-distributed but intense outgoing hux nearer the midplane—was also seen in 
the two dimensional a/M = 0.938 simulation of McKinney & Gammie (2004). 


5. Mass Outflow 
5.1. The centrifugal barrier 

We conhrm by direct examination of the complete data hies at late time that the mass 
outflow in the jets is confined to the funnel wall component. In fact within the funnel at the 
radius for which the outflow fluxes are measured, r = lOOM, the density inside the funnel has 
dropped to the floor value; matter is, of course, excluded from the funnel if it has signiflcant 
angular momentum. The mass-weighted mean angular momentum in the funnel-wall outflow 
can be computed as 

{U^)p = j dV^pU^W/ j dV^pW. (7) 

This quantity is always comparable to the angular momentum of Keplerian orbits in the 
innermost part of the disk. For example, in KDPg, {U^)p rises from ~ 1 at the base of the 
outflow (r ~ 3M) to ~ 2 at the outer boundary (see Fig. 6); by comparison, U^^rms) = Ims = 
2.1. At higher spin, {U^)p can rise as high as ~ 3, or roughly twice /ms- 
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Fig. 5.— Time-averaged Poynting flux on the inner boundary in code units times 10®. Left 
panel is KDIb, averaged over 26 snapshots equally spaced from t = 8000M to lOOOOM. The 
dashed line shows a simple fit to sin^ 9. Right panel is KDH, averaged over 27 snapshots 
from t = 7920M to lOOOOM. 




Fig. 6.— Specific angular momentum and mass flux at late time in the KDPg simulation. 
Color contours show {U^)p. Superimposed line contours are logarithm of the mass flux in 
the outflow. Peak mass flux follows very closely along the {U^)p ~ 2 contour. 
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Because the mean specific angular momentum changes little through the outfiow, it 
makes sense to describe the outfiow dynamics as taking place with respect to a single effective 
potential 


U = <! + 


9 


M 




( 8 ) 


That is, U is the energy Ut a test particle would have if orbiting with = {U^)p and 
Ur = Ug = 0 at a given location. We display the correspondence between the mass outflow 
and the shape of this effective potential in Figure 7. Two examples are presented, a/M = 0.9, 
for which — 2.5, and a/M = —0.9, for which {U^)p = 4. In both, the mass flow must 
climb in the effective potential out to r ~ lOM, but “slides” at roughly constant potential 
from there outward. Although illustrated only for these two cases, this pattern is generally 
seen in the other high-spin simulations as well. 
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Fig. 7.— Mass outflow {\/—gpU^) and effective potential contours at the same late time in 
simulation KDPg (left panel) and KDR (right panel) as in Fig. 6. Note the expanded scale 
of this hgure relative to that one. Here, the color contours are logarithm of the mass outflow 
rate, while the line contours are linear in the potential with separation 0.025. The mean 
angular momentum used to dehne the effective potential is 2.5 in the left panel and 4 in the 
right. 
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It should not be too surprising that, for the most part, the outflow moves through 
the region of least effective potential permitting access to inhnity. Thus, the shape of the 
mass outflow (at least in the inner few tens of gravitational radii) is determined by two 
factors: its specihc angular momentum (which is the angular momentum of the inner disk) 
and the gravitational potential. A corollary of the importance of the effective potential is 
that there is very little mass outflow (or mass) at polar angles closer to the axis than the 
lowest open effective potential contour. The angular momentum of the matter creates a 
centrifugal barrier that prevents matter from hlling in the funnel. Since the only source of 
mass is the accretion disk, the funnel is empty. In general, provided there is little matter of 
truly small angular momentum anywhere in the system, we expect that the mass flux in all 
such outflows should trace out a hollow cone, while the interior of the cone should have very 
low density. 

The sharp density gradient in polar angle at the inner edge of the funnel wall outflow 
poses a difficulty for numerical simulations: an extremely fine gridscale is required in the 
vicinity of the centrifugal barrier. The simulations reported here, while well-resolved in terms 
of density elsewhere, fail to describe the sharpness of this gradient. The fractional contrast 
in density from one cell to the next (moving purely in 6) on the inner edge of the funnel wall 
is typically ~ 1-10. Because of this poor resolution at the edge of the funnel, the rest-mass 
and internal energy content of the funnel interior may not be well-described quantitatively. 
Because there is negligible matter in the simulation with small enough angular momentum to 
enter the funnel interior, anything found there came either from invocation of the numerical 
density floor, or via numerical diffusion that resulted in a low ratio of angular momentum to 
matter. Of course this process of angular momentum “filtration” should exist in real systems, 
since only very low angular momentum material is allowed within the funnel. Unfortnnately, 
it is difficnlt to estimate how mnch truly low angular momentum matter may be found close 
to the center of an accretion flow—so little mass is needed in order to snbstantially change 
the content of the fnnnel interior that many physical processes, justly neglected in the rest 
of the simulation, may contribnte. For these reasons, we cannot be quantitative regarding 
the matter content in the fnnnel interior. 


5.2. Mass outflow rate 

Althongh the mass density inside the fnnnel is ill-determined, the mass density in the 
fnnnel wall, where snbstantially all of the nnbonnd mass ontflow takes place, is mnch more 
reliable. We can therefore examine how the rate of mass loss in the outflow compares to the 
accretion rate, and how this changes as a fnnction of black hole spin. 



However, there are certain cautions that attach even to this quantity. First, the ratio of 
mass outflow rate to mass accretion rate is meaningful only when averaged over long time 
intervals that exclude the initial transients. On the one hand, the accretion rate into the 
black hole is highly variable, and the timescales of these variations are often quite short, 
while on the other hand, the outflow varies more slowly because its characteristic timescale 
f'ont/c is comparatively long. In addition, there are large amplitude transients in both the 
mass inflow and outflow rates that typically occur ~ 2000M after the beginning of each 
simulation. We therefore consider the ratio between the total mass lost in in the outflow to 
the total rest-mass accreted over the time 3000-10^M. 

Second, although the density gradient within the funnel wall outflow is modest, its 
angular width in terms of cells is only about half a dozen. Thus, there remains a noticeable 
numerical uncertainty in any quantity integrated over the funnel wall. 

Third, the distinction between unbound and bound mass flow is a subtle one. Wherever 
there is signihcant mass along the outer boundary, it is moving outward. The specific energy 
{—hUt) varies smoothly with polar angle at this boundary, so that while one can dehne a 
boundary at —hUt = 1, there is no sharp dividing line in energy between the unbound and 
bound outflow near the funnel wall. 

Finally, as shown in Paper III, the mean unbound mass outflow rate increases with 
increasing radius. Beyond r ~ 20M, the scaling is slow (roughly oc but it is clearly 

still varying within our simulation volume. Consequently, we stress that the numbers given 
are the mass flows at a fiducial large radius, and not necessarily the mass flow that would 
reach inhnity. 

The mass outflow data, subject as they are to all these caveats, are given in Table 
2. Very little mass is expelled in the Schwarzschild black hole model. When a/M exceeds 
0.9, the rate begins to rise sharply. It is noteworthy that the mass loss rate relative to the 
accretion rate in the retrograde case, a/M = —0.9, is almost the same as in the prograde case 
with the same magnitude spin parameter. There is little apparent trend between a/M = 0.9 
and 0.95; for all these simulations, the mass loss rate is ~ 8% of the accretion rate. The 
scatter between these numbers is likely an estimate of (or perhaps a lower bound on) their 
uncertainty. At the highest spin, 0.99, the mass lost in the funnel wall outflow is nearly a 
quarter of the mass accreted by the black hole. 

Similarly, the average mass outflow velocity is at least partially a function of black hole 
spin. We compute an average velocity by the ratio of the time-averaged funnel wall mass flux 
to the time-averaged density in the funnel wall jet. As with many quantities, the value of the 
velocity is a function of radius. For the a/M = 0 model, the mean funnel wall jet velocity 



is w/c = 0.24 at r = lOOM. The retrograde spin model has a comparable jet velocity, while 
the velocity at r = lOOM increases io v/c = 0.3-0.4 for the high spin models. 


5.3. Outward forces 

We now embark on a search, at a qualitative level, for the forces responsible for driving 
this outflow. The most elementary evidence in this regard can be seen in Figure 8, which 
shows the azimuthally-averaged matter momentum outflow in the funnel-wall, i.e., 

V = (2/7r) j d(t) y/^pU^UrW. (9) 

We omit the factor h from the expression for the momentum density in order to suppress 
contributions from the funnel interior; h ~ 1 everywhere in the funnel wall. As this hgure 
shows, the matter momentum outflow in the funnel wall increases out to r ~ 20M. Time- 
dependent fluctuations cause the momentum outflow to vary with increasing radius at larger 
distance, but without any consistent trend. Moreover, the track of the momentum outflow 
follows a path along which the effective potential hardly changes. Figure 6 also provides 
evidence pointing toward the same conclusion. The hgure illustrates that the mean specihc 
angular momentum in the funnel-wall hardly changes with radius. Thus, the effect of any 
torques must be quite minor. 

These results indicate that while centrifugal effects help guide the outhow, they do not 
accelerate it: other forces must be responsible for its radial acceleration. To ascertain just 
what these forces are, we look at the problem from a different point of view. In Paper III, 
we demonstrated that mass in the outhow is acquired rapidly in the range between a few 
and ~ 20M, but continues to grow slowly with increasing radius beyond that point. In order 
for mass to be introduced to the outhow over such a broad range of radii, the acceleration 
region must be distributed along the length of the funnel wall. Just such a situation can be 
seen in Figure 9. In that hgure, we show a Newtonian approximation to the Mach number 
relative to the magnetosonic speed, i.e., where = [||&|P J- (5/3)p]/(| |6|p J- hp). 

Regions outside the unbound outhow are suppressed in this hgure. From r ~ lOM to the 
outside, the transonic contour follows the outside edge of the funnel-wall. That is, fresh 
matter is both injected and accelerated along the length of the jet. In much of the funnel 
interior, the how is sub-magnetosonic because ||&|P/p is so large that Cms — 1; at larger radii 
in the funnel interior, R’’ —> 1, so that, in this sense, the how approaches Mach number 
unity. 
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Fig. 8.— Momentum outflow V and effective potential contours at late time in simulation 
KDH {a/M = 0.95). The color contours are linear in P, while the line contours are linear in 
the potential for = 3.0 (the most appropriate for this simulation) with separation 

0.05. 
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Fig. 9.— Azimuthal average of the approximate magnetosonic Mach number in matter 
belonging to the outflow at a late time in the a/M = 0.5 simulation. 
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Pressure forces account for this distributed acceleration. Figure 10 shows the total 
pressure (more precisely, ^J—g{\\h\\^/2 +p)) compared to the radial matter momentum flux 
in the outflow {y/—gphU^Ur). The pressure contours are oblique to the momentum flux 
contours, angled at such a way as to push matter into the outflow as well as to accelerate 
it radially. Matter is prevented from moving deeper into the funnel interior by the sharp 
rise in the effective potential (cf. Fig. 8). In the coronal region generally, and near the 
funnel wall in particular, the magnetic and gas pressures are comparable to one another. 
Interestingly, the contours of both the gas pressure and the magnetic pressure considered 
alone are signihcantly more irregular than the contours of their sum. 


5.4. Collimation 

We have already demonstrated that the inner edge of the matter outflow is determined 
by the centrifugal barrier in the inner few tens of gravitational radii, and that its flow 
beyond that distance is essentially ballistic within the nearly-flat effective potential. The 
electromagnetic helds in the funnel interior expand to £11 the hollow region, so the matter’s 
dynamics determine the edge of the electromagnetically-dominated funnel interior. The 
outer edge of the matter outflow (funnel wall jet) is less sharply marked. The specific energy 
{—hUt) is a smooth function of polar angle at large radius, and there is no sharp division 
between bound and unbound regions. Similarly, the pressure is likewise smoothly varying 
(cf. Paper III). In this sense, we consider the outflow to be collimated by the coronal 
pressure. For this reason, future work with more realistic coronal thermodynamics could 
reveal interesting changes in outflow shape. However, the smoothness of the border between 
outflow and bound corona also means that it is ill-marked. Indeed, the rise of the mass 
outflow rate with radius may be an indication that entrainment plays a significant role in 
mass entry at large radius, further blurring this edge. 

To quantify the collimation, we plot in Figure 11 the opening angle of the jet obtained 
by locating the angular boundary between the bound and unbound flow as a function of 
radius. As the figure shows, the higher the spin, the more narrow the jet initially, although 
all of the jets collimate out to around r = 40-50M. Beyond that point, the opening angle 
remains constant or increases very slightly toward the outer boundary of the simulation 
domain. The spin-dependence of the opening angle at small radii is most likely due to the 
fact that the stable disk extends farther inward with increasing spin. Where there is stable 
disk in the equatorial plane, it is easier to support a high-pressure corona. 
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Fig. 10.— Azimuthal average of the total pressure (gas plus magnetic) in color contours 
at the same time in the a/M = 0.5 simulation shown in Fig. 9. Radial momentum flux in 
the outflow is shown in white contours. Both sets of contours are linear scales. The highest 
contour in momentum flux is roughly twice the highest pressure shown in the color scale. 





Fig. 11.— Opening angle of the jet in radians, as defined by the angular location of the 
boundary between bound and unbound outflow. Four models are shown at t = 10, DOOM: 
the solid line is a/M = —0.9, the dotted line is a/M = 0.5, the dashed line is a/M = 0.9, 
and the dot-dashed line is a/M = 0.99. 
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6. Summary and Discussion 

In this paper we have examined the unbound axial outflows (jets) that develop within 
GRMHD accretion simulations. What is noteworthy in these models is that the jets develop 
as a natural result of the accretion process, rather than as a direct consequence of conditions 
imposed on the simulation. The jets have two major components: a matter-dominated 
outflow that moves at a modest velocity {vjc ~ 0.3) along the centrifugal barrier surrounding 
an evacuated axial funnel, and a highly-relativistic Poynting flux dominated jet within the 
funnel. The funnel wall jet is accelerated and collimated by magnetic and gas pressure 
forces in the inner torus and the surrounding corona. The Poynting flux jet results from 
the formation of a large scale radial magnetic held within the funnel. This held is spun by 
the rotating spacetime of the black hole, hence the energy for this jet component comes not 
directly from accretion but from the black hole’s rotation, although the held mediating this 
energy extraction is a product of accretion. The energies in both components of the jet are 
a function of black hole spin; greater spins yield greater energy. The total jet power can be 
a significant fraction of the accretion energy estimated by the product of the mass accretion 
rate and the specific binding energy of the last stable orbit. 

We next consider the broader implications of these results within the context of several 
of the standard theoretical categories for jets, namely jet acceleration, power, collimation, 
and the Blandford-Znajek mechanism. 


6.1. Accelerating forces 

It has been expected for many years (since Blandford & Znajek 1977), that within the 
axial funnel the spacetime rotation enforced by spinning black holes can drive a significant 
electromagnetic outflow. Our results are clearly in line with that expectation. 

The outflow-generation process within the accretion flow itself that has for many years 
received the greatest attention has been magnetocentrifugal acceleration, such as envisioned 
in the model of Blandford & Payne (1981), but this mechanism does not appear to be 
important here. 

On the other hand, the main driving force we identify for the funnel wall jet—a high- 
pressure corona squeezing material against an inner centrifugal wall—has received very little 
attention in the past. In consequence of this compression, the pressure at the base of the 
matter outflow is kept high enough to create an outward radial pressure gradient that pro¬ 
vides the momentum for the outflow. It is also worth commenting that these forces act over 
a substantial radial range, from a few gravitational radii out to many tens. 



6.2. Outflow strength 


To measure the relative power in the unbound outflow, we have compared the jet energy 
flux with the mass accretion rate into the hole over the last 70% of the simulation. The 
simulation has no radiative losses, so we cannot say what the radiative efficiency of the 
accretion disk would be. The Novikov-Thorne efficiency provides a convenient standard of 
comparison. This nominal efficiency corresponds to the binding energy of the innermost 
stable circular orbit, although the actual radiative efficiency could well be different because 
of, for example, magnetic torques in the plunging region (Krolik 1999; Gammie 1999), or 
heat advection into the hole. Judged by this measure, the energy carried off in the outflow is 
considerable. In two cases (a/M = —0.9 and 0.99), the electromagnetic power at r = lOOM 
is, by itself, comparable to the Novikov-Thorne prediction of radiated power. In the other 
high spin cases, it is about a quarter of the Novikov-Thorne number. One should remember, 
however, that energy can be (and often is) exchanged between the fields and the matter, 
so these measures of the electromagnetic power in the outflow are quantitatively accurate 
only at the fiducial location where they are measured. Because the rate of energy loss in the 
matter outflow is usually a few times as great as the electromagnetic contribution, the total 
is typically comparable to the classical prediction of radiated power. 

Another way to gain a feel for the magnitude of the electromagnetic power in the 
outflow is via dimensional analysis. The Poynting flux at the inner boundary is oc 
(cf. eqn. 6); but what determines the magnetic held strength? Indeed, this has been the 
key question examined by such works as Ghosh & Abramowicz (1997) and Livio, Ogilvie, 
& Pringle (1999). These authors argued that the strength of the held responsible for a 
Blandford-Znajek process must be set by the fields in the accretion flow, and this may limit 
the total jet power. 

At radii in the disk well beyond the marginally stable region, angular momentum 
conservation ensures that oc MQ/h if magnetic torques are the primary agent of 

accretion, while the competition between orbital shear and other stresses tends to make 
(||6|p) oc In these simulations, in which h/r is common to all, but the black hole 

spin varies, if we wish to analyze the spin-dependence of the Poynting flux near the event 
horizon, it makes sense to examine the ratio of the magnetic held intensity to the mass ac¬ 
cretion rate at that location. We And that the time- and surface area-averaged ratio ||6|p/M 
grows rapidly with increasing black hole spin, very nearly oc (1 — \a/M\)~^'^ (Fig. 12). Some 
of this growth is attributable to mere compression, as the event horizon shrinks as the spin 
rate rises. However, the dependence is much stronger than anything simple geometry could 
explain. It appears that the strong frame-dragging associated with extreme Kerr spacetimes 
leads to substantial held amplification as well. The rapid growth in electromagnetic power 



output with faster black hole rotation is predominantly due to this effect. 
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6.3. Jet shape 

One of the key issues in jet formation and propagation is the degree to which the jets 
are collimated. Of course, jets are observed on far larger scales than we are modeling in 
these simulations, but the degree of collimation seen here is nevertheless of some interest. 
Generally speaking, whatever collimation there is occurs in the inner regions of the flow, 
where the inner disk and hot corona provide confinement. Relatively close to the black hole 
(r < lOM), the boundary between the magnetically-dominated Poynting flux jet and the 
funnel wall outflow follows the shape of the centrifugal barrier. At larger distances, however, 
the matter flows almost straight outward, so that the funnel is nearly conical, while the 
contours of constant effective potential gradually move toward smaller polar angle (Fig. 8). 

The precise shape and collimation of the jet is somewhat uncertain because the outer 
boundary of the funnel-wall, matter-dominated funnel wall outflow is somewhat indistinct. 
We see a smooth transition as a function of polar angle between mildly relativistic outflowing 
unbound matter and slightly slower, but bound, coronal matter. On the other hand, the 
boundary between the low-density magnetically-dominated funnel interior and the higher- 
density funnel-wall outflow is sharp and clear. It is possible that magnetic forces might 
provide additional collimation on larger scales, but this remains to be seen in future simula¬ 
tions. 


6.4. Comparison to other models and the “Blandford-Znajek” mechanism 

In the astrophysical literature, the general idea that energy can be conveyed from rotat¬ 
ing black holes to the outside world via magnetic connections is often dubbed the “Blandford- 
Znajek mechanism” because they were the hrst to discuss a specihc version of this process. 
In our simulations, when all is said and done, the funnel interior is dominated by a magnetic 
held that would be purely radial but for the rotation created by its connection to the black 
hole’s ergosphere. Thus, any toroidal held in the funnel, and consequently any Poynting 
hux, is due entirely to general relativistic dynamics peculiar to spinning black holes. In 
other words, the outward Poynting hux must be due to a process that could ht under the 
general rubric of “Blandford-Znajek”. 

On the other hand, we also hud specihc diherences from the classical Blandford-Znajek 
model, specihcally the particular idealized picture presented in their original paper (Bland- 
ford & Znajek 1977). For example, although the funnel region is magnetically-dominated, it 
is not in general in a state of force-free equilibrium. Indeed, the very large huctuations that 
continually occur in the outhow show that it is never in any state of equilibrium, force-free 



- 35 - 

or otherwise. In addition, as originally shown in Phinney (1983), the collapse of all wave 
speeds to c in the exact force-free limit means that to understand the dynamics properly 
requires retention of a hnite inertia. 

Moreover, whereas in the classical Blandford-Znajek model there was zero accretion (the 
only function of nearby orbiting matter is to support the currents sustaining the magnetic 
held), in our picture the electromagnetic luminosity in the outhow takes place entirely as a 
byproduct of accretion. In the classical model the magnetic held was a relic of plasma inhow 
that took place in the distant past; because there was no continuing accretion, the energy 
output could be attributed unambiguously to the black hole’s rotation. In our simulation, 
the magnetic held is a direct result of on-going accretion. While the proximate energy source 
for the jet is the black hole’s rotation, accretion replenishes both the black hole’s mass and 
angular momentum, and whether the rotational energy of the black hole actually decreases 
depends on details such as the radiative efficiency and, possibly, the accretion of photons 
(Krolik 2001; cf. Thorne 1974). In our a/M = 0.99 simulation, the net angular momentum 
accreted is so small that, for reasonable values of the radiative efficiency, the system is very 
near the edge of an actual decrease in rotational kinetic energy. Even if the black hole’s 
rotational energy does not decrease, there can nevertheless be a significant increase in the 
efficiency of accretion. We hud that for rapidly rotating black holes the total power in the 
Poynting flux jet can be comparable to the nominal accretion power. 

A number of other numerical studies of the generic Blandford-Znajek mechanism have 
been carried out recently. One type of simulation has been devoted to examining the idealized 
problem of a rotating black hole embedded in thin plasma containing a specihed external 
magnetic held (Koide 2003; Semenov et al. 2004; Komissarov 2005). These simulations 
focus on points of basic principle, such as the extent to which accretion of either matter 
or electromagnetic helds with negative energy-at-in£nity plays a role. Komissarov (2005) 
particularly notes that plasma effects within the ergosphere were significant early on in his 
simulation while the initial plasma is accreting into the hole, but that the long-term steady 
state mechanism resembles closely the classic Blandford-Znajek model. 

The simulations of McKinney & Gammie (2004) were more similar in their focus to ours, 
as they too were principally concerned with the connection between continuing accretion 
and the electromagnetic power output of black holes. For this reason, their (axi-symmetric) 
simulations employed initial conditions resembling ours (isolated initial gas torus containing 
held loops). They found good agreement with several aspects of the classical Blandford- 
Znajek model for slowly rotating holes, but growing differences as the spin of the hole 
increases and the held becomes less monopolar. That this should be so is not too surprising, 
as the specific predictions of the Blandford-Znajek model were developed on the basis of 



a perturbation expansion in a/M. In regard to the global energy and angnlar momentnm 
bndgets, McKinney & Gammmie’s results qualitatively resemble ours. They also saw that 
as the hole loses energy and angnlar momentnm electromagnetically, its total mass and 
angular momentum are replenished by accretion. Likewise, they, too, fonnd that the overall 
efficiency (energy released to inhnity per rest-mass accreted) can be significantly enhanced by 
the energy carried off in the jet. McKinney (2005) expanded npon these resnlts and proposed 
an analytic fitting fnnction for the electromagnetic efficiency of the jet as a fnnction of black 
hole rotation rate, snggesting that this efficiency scales oc where VLh is the rotation rate 
of the black hole. If the Poynting flux were due simply to rotation of a magnetic dipole 
whose magnitude is independent of black hole rotation bnt oc the efficiency would, 

of course, increase oc The slightly faster dependence McKinney (2005) fonnd snggests 
either that there is a mixtnre of a higher multipole component or that the magnitude of the 
magnetic dipole increases slowly with increasing black hole spin. 

The choice of this fnnctional form was motivated by the hope that electromagnetic 
radiation in this case resembled classical dipole radiation, in which the Inminosity is oc 
Whether this applies in the present context depends, of course, on whether dipole field 
geometry is a reasonable approximation and on whether the mean field strength at the 
horizon prodnced by a fixed accretion rate is independent of spin. 

Unfortnnately, for a variety of technical reasons it is difficult to compare quantitatively 
the results of different nnmerical simnlations. Withont the inclusion of additional physics 
in the basic equations, there is no absolute density scale in the gas. Althongh efficiencies 
are well-defined (ratios of energies to masses within the simnlation), one cannot define an 
absolnte accretion rate. Even the comparison of efficiencies from one type of simulation to 
another poses problems. For example, two dimensional simnlations tend to over-emphasize 
axisymmetric MRI modes. As a resnlt, they tend to have large accretion rates early on, 
followed by a rapid decline as the magnetic field and accompanying tnrbnlence die ont, 
as expected from the anti-dynamo theorem. Thus, when compnting the efficiency in a 
two dimensional simulation simulation, one must be careful about choosing an appropriate 
interval over which to average the accretion rate. A second problem is that, as discussed 
previously, energy can be exchanged between the electromagnetic and matter components of 
the ontfiow, so the pnrely electromagnetic power is a fnnction of radius. To cross-compare 
different simnlations, it is therefore important to use the same fiducial distance for evalnating 
the electromagnetic power. Lastly, jet properties may depend on the vertical thickness of 
the accretion disk, a parameter determined by the initial temperatnre of the disk and the 
treatment of the energy eqnation in the simnlation. The gross accretion rate in general scales 
oc bnt we do not know how the magnetic field in the vicinity of the black hole scales 
with disk temperature. Further simulational work will be required to determine snch effects. 
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although because the strength of the axial funnel field is determined by the accretion flow, 
it is a reasonable expectation that the jet power might scale with the accretion rate (Ghosh 
& Abramowicz 1997). 

Keeping these caveats in mind, we can compare our electromagnetic efficiency in the 
jet with the scaling relation proposed by McKinney (2005) by modifying his relation to 
use the rotation rate of our inner grid radius instead of the rotation rate of the horizon. 
The quantitative match between our data and McKinney’s formula is not particularly good, 
although qualitatively we, too, see a strong increase in efficiency as black hole spin increases. 
Our results are better described by the simple formula 

~ 0.002/(1 - |a/M|), (10) 

a form that describes both the high and lower spin models. Of course, without any clear 
theory to predict the strength of the magnetic field on the event horizon relative to the mass 
accretion rate as a function of black hole rotation rate, one cannot expect to find a general 
formula for total jet power. 

Finally, we would like to close this paper by pointing out that, although the Blandford- 
Znajek mechanism was originally thought of as a scheme for powering jets, outflows are not 
the only place where electromagnetic power from the black hole can go. In this paper we 
have concentrated on jet properties. However, preliminary analysis of the accretion flow in 
these simulations indicates that the Poynting flux and electromagnetic angular momentum 
flux from the black hole into the disk can be comparable to those within the outflow. These 
fluxes and their consequences will be examined in subsequent work. 


This work was supported by NSF grant PHY-0205155 and NASA grant NNG04GK77G 
(JFH), and by NSF grants AST-0205806 and AST-0313031 (JHK). We acknowledge Jean- 
Pierre De Villiers for continuing collaboration in the development of the algorithms used 
in the GRMHD code. The simulations were carried out on the DataStar system at SDSG. 
A portion of this work was carried out while the authors were at the Kavli Institute for 
Theoretical Physics, supported by NSF grant PHY99-07949. 

REFERENCES 

Balbus, S. A. & Hawley, J. F. 1998, Rev. Mod. Phys., 70, 1 
Blandford, R. D. & Znajek, R. 1977, MNRAS, 179, 433 
Blandford, R. D. & Payne, D. 1981, MNRAS, 199, 833 



De Villiers, J.-P. & Hawley, J. F. 2003, ApJ, 589, 458 

De Villiers, J.-P., Hawley, J. F., & Krolik, J. H. 2003, ApJ, 599, 1238 (Paper I) 

De Villiers, J.-P., Hawley, J. F., Krolik, J. H. & Hirose, S. 2005, ApJ, 620, 878 (Paper HI) 
De Villiers, J.-P., Staff, J., & Ouyed, R., astro-ph 0502225 
Evans, C. R., & Hawley, J. F., ApJ, 332, 659 
Gammie, C. F. 1999, ApJ Lett 522, L57 

Gammie, G. F., McKinney, J. G., & Toth, G. 2003, ApJ, 589, 444 
Ghosh, P. & Abramowicz, M. A. 1997, MNRAS, 292, 887 

Hirose, S., Krolik, J. H., De Villiers, J. P., & Hawley, J. F. 2004, ApJ, 606, 1083 (Paper H) 

Kato, Y., Mineshige, S., & Shibata, K. 2004, ApJ, 605, 307 

Koide, S. 2003, Phys Rev D, 67, 104010 

Komissarov, S. S. 2004, MNRAS, 350, 1431 

Komissarov, S. S. 2005, MNRAS, 359, 801 

Krolik, J. H. 1999, ApJ Lett 515, L73 

Krolik, J. H. 2001, in Explosive Phenomena in Astrophysical Compact Objects, 1. Yi and M. 
Rho, eds. (AIP: New York), p. 10 

Krolik, J.H., Hawley, J.F. & Hirose, S. 2005, ApJ 622, 1008 (Paper IV) 

Livio, M., Ogilvie, G. L, & Pringle, J. E. 1999, ApJ, 512, 100 

Lynden-Bell, D. 2003, MNRAS 341, 1360 

McKinney, J. G. 2005, ApJ, 630, L5 

McKinney, J. G., & Gammie, G. F. 2004, ApJ 611, 977 

Novikov, I.D. & Thorne, K.S. 1973, in Black Holes: Les Astres Occlns, eds. G. de Witt & B. 
de Witt (New York: Gordon & Breach), 344 

Phinney, E.S. 2003, nnpnblished Gambridge University Ph.D. thesis 

Pnnsly, B., & Goroniti, F. V. 1990, ApJ, 354, 583 



Semenov, V., Dyadechkin, S., & Punsly, B. 2004, Science, 305, 
Shibata, K. & Uchida, Y. 1985, PASJ, 37, 31 
Stone, J. M., & Norman, M. L. 1992a, ApJS, 80, 753 
Stone, J. M., & Norman, M. L. 1992b, ApJS, 80, 791 
Thorne, K. S. 1974, ApJ, 191, 507 


This preprint was prepared with the AAS IA]J]X macros v5.2. 




Table 1. Simulation Parameters 


Name 

a/M 

r+/M 

nn/M 

hPmax 

Cc;(nn) 

rms/M 

U^{rms)lM 

KDR 

-0.9 

1.44 

1.503 

10. 

-0.287 

8.72 

4.17 

KDOc 

0.0 

2.00 

2.104 

10. 

0.0 

6.00 

3.46 

KDIb 

0.5 

1.86 

1.904 

10. 

0.127 

4.23 

2.90 

KDPg 

0.9 

1.44 

1.503 

10. 

0.289 

2.32 

2.10 

KDG 

0.93 

1.37 

1.458 

6. 

0.306 

2.10 

1.98 

KDH 

0.95 

1.31 

1.403 

6. 

0.326 

1.94 

1.89 

KDJa 

0.99 

1.14 

1.203 

6. 

0.406 

1.45 

1.57 


Note. — Here a/M is the spin parameter of the black hole, r+ is the 
horizon radius, rjn is the innermost radius in the computational grid, 
hPmax is the largest Lorentz factor permitted, u}{r-) is the frame-dragging 
rate at that radius in the equatorial plane (i.e., the rotation rate of a 
zero angular-momentum observer), r^s is the radial coordinate of the 
innermost stable circular orbit in the equatorial plane, and is 

the specific angular momentum of that last stable orbit. 






Table 2. Jet Mass Flux 


Model 

a/M 

Mtot 

Afhole 

Mjet 

Mjet/Mhoie 

KDR 

-0.90 

195.6 

22.1 

0.77 

0.035 

KDOc 

0.00 

156.3 

22.2 

0.11 

0.005 

KDIb 

0.50 

234.8 

20.6 

0.76 

0.037 

KDPg 

0.90 

290.5 

19.2 

1.71 

0.089 

KDG 

0.93 

319.0 

21.8 

1.72 

0.079 

KDH 

0.95 

339.2 

20.3 

1.72 

0.085 

KDJ 

0.99 

382.5 

15.5 

3.66 

0.237 


Note. — Mtot is the total mass in the initial hy¬ 
drostatic torus; Mhoie is the total rest-mass passing 
through the inner boundary in the time interval 
3000-10000M; Mjet is the total unbound rest-mass 
leaving the simulation region over the same inter¬ 
val. All masses are given in the (arbitrary) code 
mass-units. 








Table 3. Angular Momentum Flux 


Model 

a/M 

r m 

hole 

r em 

^ hole 

r m 

^ jet 

r em 
^ jet 

KDR 

-0.90 

-82.7 

-4.8 

3.15 

-2.65 

KDOc 

0.00 

-72.4 

1.4 

0.31 

0.14 

KDIb 

0.50 

-57.0 

2.6 

1.40 

2.10 

KDPg 

0.90 

-42.9 

7.5 

2.26 

6.96 

KDG 

0.93 

-46.0 

11.7 

5.44 

6.36 

KDH 

0.95 

-42.2 

16.8 

5.68 

8.88 

KDJ 

0.99 

-28.7 

26.4 

14.3 

12.8 


Note. — L'^oXe is total angular momentum 
carried by the mass passing through the inner 
boundary in the time interval 3000-10000M; 
is the angular momentum carried by the electro¬ 
magnetic held. The values marked “jet” are the 
total angular momenta leaving the simulation re¬ 
gion over the same interval in the unbound jet 
outhow. 







Table 4. Energy Flux 


Model 

a/M 

jT'm 

'^hole 

rpem 

■^hole 

JT'm 

-^jet 

Vm 

jrem 

-^jet 

Vem 

hNT 

KDR 

-0.90 

-20.8 

0.41 

2.71 

0.088 

0.509 

0.023 

0.039 

KDOc 

0.00 

-20.6 

-0.02 

0.16 

0.0022 

0.007 

3.1 X 10“^ 

0.057 

KDIb 

0.50 

-21.7 

0.10 

2.06 

0.063 

0.129 

0.0063 

0.081 

KDPg 

0.90 

-21.8 

0.88 

5.89 

0.22 

0.892 

0.046 

0.155 

KDG 

0.93 

-19.8 

1.37 

3.13 

0.065 

0.824 

0.038 

0.173 

KDH 

0.95 

-16.9 

2.79 

4.26 

0.13 

1.46 

0.072 

0.190 

KDJ 

0.99 

-9.9 

6.86 

9.94 

0.41 

3.28 

0.21 

0.264 


Note. — Negative signs for Ehoie signify energy flows into the black hole; 
positive signs for that quantity indicate energy flowing into the compu¬ 
tational domain. For dehnitions of the several efficiencies denoted by rj 
with different subscripts, see the text. 







